*************************************************************************************************
*																								*
*						Careers - Rambachan and Roth Honest DiD test (JPEMicro R1)				*
*						name file: rambachan_roth.do											*
*																								*
*************************************************************************************************

** NOTES


* RM (Relative magnitude)
/* "In empirical applications, researchers may be willing to assume that the confounding factors which produce non-parallel trends in the posttreatment periods are not too much larger in magnitude than the confounding factors in the
pre-treatment periods. Such a restriction can be
formalized in our framework by imposing that δ ∈ DELTA RM (M¯ ) for M¯ ≥ 0, where
DELTA RM (M¯ ) = {δ : ∀ t ≥ 0, |δt+1 − δt| ≤ M¯ · max
s<0
|δs+1 − δs|}.  The choice DELTA RM (M¯ ) may be reasonable if the researcher
suspects that possible violations of parallel trends are driven by confounding economic shocks
that are of a similar magnitude to confounding economics shocks in the pre-period. When
the number of pre-treatment and post-treatment periods is similar, a natural benchmark may
be M¯ = 1, which bounds the worst-case post-treatment difference in trends by the equivalent
maximum in the pre-treatment period."
*/


* SD (Smoothness restrictions) - THIS SEEMS MORE APPROPRIATE IN OUR CONTEXT
/*
In other empirical settings, researchers may be worried
about confounding from secular trends (e.g. long-run changes in labour supply) that they suspect
evolve smoothly over time. . In such settings, it is common for empirical researchers to control
for a linear group-specific time trend. A natural relaxation is therefore
to impose only that the differential trends evolve smoothly over time by bounding the extent to
which its slope may change across consecutive periods. Such a restriction can be formalized in
our framework by imposing that δ ∈ DELTA SD(M) for M ≥ 0, where
DELTA SD(M) := {δ : |(δt+1 − δt) − (δt − δt−1)| ≤ M, ∀ t}.  [i.e., M is the order of the slope function]

The parameter M ≥ 0 governs the amount by which the slope of δ can change between consecutive periods, and thus bounds the discrete analogue of the second derivative. In the special case where
M = 0, DELTASD(0) requires that the difference in trends be exactly linear, which corresponds with
the assumption underlying the parametric linear specification common in applied work.
*/

* BOTH

/*
When the setting lends itself to RD but we can't figure out M, we can combine the two and have SD where we select M to be the Maximum slope change in the pre-period.
*/

* Install programs -- *
ssc install honestdid, replace
* -- *

clear all
set more off
set matsize 11000
set maxvar 32000
cd "~/Dropbox/Research/entrepreneurs/dofiles/dofiles_career/revision_JPEMicro_R1/replication"
global out = "~/Dropbox/Research/entrepreneurs/draft bipolar/figures/JPEMicro_2024"
global graph = "xlabel(,grid glp(dot) glc(gs10) tlc(gs10)) ylabel(,grid glp(dot) glc(gs10) tlc(gs10)) plotregion(lp(blank))"

* Yale colors
global yaleblue = "0 53 107"
global ylb = "40 109 192"
global yo = "189 83 25"
global yg = "94 112 45"
global ylight = "217 233 242"
global graph = "xlabel(,grid glp(dot) glc(gs10) tlc(gs10)) ylabel(,grid glp(dot) glc(gs10) tlc(gs10)) plotregion(lp(blank))"
global graphbar = "ylabel(,grid glp(dot) glc(gs10) tlc(gs10)) plotregion(lp(blank))"



*--------> NO earnings, triple difference
 

*--------> 	Load vectors and define matrix names
import excel using "output/rambachan_roth/bd_matrix_triple_noearn.xlsx", sheet(beta) clear

mkmat A-B D-N, mat(A)
matrix define beta = A[1,1..13]
matrix list beta

import excel using "output/rambachan_roth/bd_matrix_triple_noearn.xlsx", sheet(V) clear

drop if _n > 14
mkmat A-O, mat(V)
matrix define V1 = V[1..2,1..2]
matrix define V2 = V[4..14,1..2]
matrix define V3 = V[1..2,4..14]
matrix define V4 = V[4..14,4..14]
matrix define V = (V1, V3 \ V2, V4)
matrix list V


*--------> 	Use command RD + SD - 1960
matrix l_vec =  0 \ 0 \ 1 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0  // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("No earnings, k = 1960") ylabel(-0.8(0.1)0.2) xsize(4)




*--------> 	Use command RD + SD - 1966
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 1 \ 0 \ 0 \ 0 \ 0 \ 0 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("No earnings, k = 1966") ylabel(-0.8(0.1)0.2) xsize(4)




*--------> 	Use command RD + SD - 1970
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 1 \ 0 \ 0 \ 0 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("No earnings, k = 1970") ylabel(-0.8(0.1)0.2) xsize(4)



*--------> 	Use command RD + SD - 1976
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 1 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("No earnings, k = 1976") ylabel(-0.8(0.1)0.2) xsize(4)



*--------> ln earnings, triple difference
 

*--------> 	Load vectors and define matrix names

import excel using "output/rambachan_roth/bd_matrix_triple_lnearn.xlsx", sheet(beta) clear

mkmat A-B D-N, mat(A)
matrix define beta = A[1,1..13]
matrix list beta

import excel using "output/rambachan_roth/bd_matrix_triple_lnearn.xlsx", sheet(V) clear

drop if _n > 14
mkmat A-O, mat(V)
matrix define V1 = V[1..2,1..2]
matrix define V2 = V[4..14,1..2]
matrix define V3 = V[1..2,4..14]
matrix define V4 = V[4..14,4..14]
matrix define V = (V1, V3 \ V2, V4)
matrix list V


*--------> 	Use command RD + SD - 1960
matrix l_vec =  0 \ 0 \ 1 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0  // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("Ln earnings, k = 1960")  xsize(4)




*--------> 	Use command RD + SD - 1966
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 1 \ 0 \ 0 \ 0 \ 0 \ 0 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("Ln earnings, k = 1966")  xsize(4)

* ylabel(-0.8(0.1)0.2)


*--------> 	Use command RD + SD - 1970
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 1 \ 0 \ 0 \ 0 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("Ln earnings, k = 1970")  xsize(4)



*--------> 	Use command RD + SD - 1976
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 1 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("Ln earnings, k = 1976") xsize(4)



*--------> No earnings, diff in difference
 

*--------> 	Load vectors and define matrix names
import excel using "output/rambachan_roth/bd_matrix_did_noearn.xlsx", sheet(beta) clear

mkmat A-B D-N, mat(A)
matrix define beta = A[1,1..13]
matrix list beta

import excel using "output/rambachan_roth/bd_matrix_did_noearn.xlsx", sheet(V) clear

drop if _n > 14
mkmat A-O, mat(V)
matrix define V1 = V[1..2,1..2]
matrix define V2 = V[4..14,1..2]
matrix define V3 = V[1..2,4..14]
matrix define V4 = V[4..14,4..14]
matrix define V = (V1, V3 \ V2, V4)
matrix list V


*--------> 	Use command RD + SD - 1960
matrix l_vec =  0 \ 0 \ 1 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0  // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("No earnings, k = 1960") ylabel(-0.8(0.1)0.2) xsize(4)




*--------> 	Use command RD + SD - 1966
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 1 \ 0 \ 0 \ 0 \ 0 \ 0 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("No earnings, k = 1966") ylabel(-0.8(0.1)0.2) xsize(4)




*--------> 	Use command RD + SD - 1970
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 1 \ 0 \ 0 \ 0 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("No earnings, k = 1970") xsize(4)



*--------> 	Use command RD + SD - 1976
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 1 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("No earnings, k = 1976") xsize(4)




*--------> ln earnings, diff in diff
 


*--------> 	Load vectors and define matrix names

import excel using "output/rambachan_roth/bd_matrix_did_lnearn.xlsx", sheet(beta) clear

mkmat A-B D-N, mat(A)
matrix define beta = A[1,1..13]
matrix list beta

import excel using "output/rambachan_roth/bd_matrix_did_lnearn.xlsx", sheet(V) clear

drop if _n > 14
mkmat A-O, mat(V)
matrix define V1 = V[1..2,1..2]
matrix define V2 = V[4..14,1..2]
matrix define V3 = V[1..2,4..14]
matrix define V4 = V[4..14,4..14]
matrix define V = (V1, V3 \ V2, V4)
matrix list V


*--------> 	Use command RD + SD - 1960
matrix l_vec =  0 \ 0 \ 1 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0  // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("Ln earnings, k = 1960")  xsize(4)




*--------> 	Use command RD + SD - 1966
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 1 \ 0 \ 0 \ 0 \ 0 \ 0 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("Ln earnings, k = 1966")  xsize(4)

* ylabel(-0.8(0.1)0.2)


*--------> 	Use command RD + SD - 1970
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 1 \ 0 \ 0 \ 0 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("Ln earnings, k = 1970")  xsize(4)



*--------> 	Use command RD + SD - 1976
matrix l_vec =  0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 0 \ 1 // size is nr of post-event periods

* Calculate maximum change in slope
local a = beta[1,1]
local b = beta[1,2]
local c = 0

disp `b' - `a'

local M = abs((`b' - `a'))
disp `M'
local M1 = `M' / 10
local M2 = `M1' * 2
local M3 = `M1' * 3
local M4 = `M1' * 4

honestdid, b(beta) vcov(V) numpre(2) delta(sd) mvec(0 `M1' `M2' `M3' `M4' `M') l_vec(l_vec) coefplot alpha(0.1) xtitle(M) xlabel(1 "Original" 2 "0" 3 "0.2" 4 "0.4" 5 "0.6" 6 "0.8" 7 "1") title("Ln earnings, k = 1976") xsize(4)

